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ABSTRACT 

The three-dimensional Monte Carlo photoionization code Mocassin has been applied to con- 
struct a realistic model of the planetary nebula NGC 3918. Three different geometric models 
were tried, the first being the biconical density distribution already used by Cle22 et al. ( 1987p. 
In this model the nebula is approximated by a biconical structure of enhanced density, em- 
t:;;^- , bedded in a lower density spherical region. Spindle-like density distributions were used for 

■ the other two models (models A and B). Model A used a mass distribution slightly modi- 
C3 I fied from one of Mellema's (1996) hydrodynamical models that had akeady been adopted by 

■ ^orradi et al. ( 1999) for their observational analysis of NGC 3918. Our spindle-like model B 
' instead used an analytical expression to describe the shape of the inner shell of this object as 

consisting of an ellipsoid embedded in a sphere. 

The effects of the interaction of the diffuse fields coming from two adjacent regions 
of different densities were investigated. These are found to be non-negligible, even for the 
relatively uncomplicated case of a biconical geometry. We found that the ionization structure 
of low ionization species near the boundaries is particularly affected. 
^1 It is found that all three models provided acceptable matches to the integrated nebular 

optical and ultraviolet spectrum. Large discrepancies were found between all of the model 
predictions of infrared fine-structure line fluxes and ISO SWS measurements. This was found 
k>( \ to be largely due to an offset of wl4 arcsec from the centre of the nebula that affected all of 

. the ISO observations of NGC 3918. 

' For each model, we also produced projected emission-line maps and position-velocity 

diagrams from synthetic long-slit spectra, which could be compared to recent HST images 
and ground-based long-slit echelle spectra. This comparison showed that spindle-like model B 
provided the best match to the observations. Although the integrated emission line spectrum 
of NGC 3918 can be reproduced by all three of the three-dimensional models investigated 
in this work, the capability of creating projected emission-line maps and position-velocity 
diagrams from synthetic long-slit spectra was found to be crucial in allowing us to constrain 
the structure of this object. 

Key words: planetary nebulae: individual: NGC 3918 - ISM: abundances 



1 INTRODUCTION 

The southern planetary nebula NGC 3918 (PN G294.6+04.7) is 
a very well known and widely studied object. A detailed study, 
based on UV, op tical and radio observations, was presented by 
IClegg et alj|l987l from now on C87). The morphological and kine- 
matical information available at the time of their work, was, how- 
ever, very limited, and, based on this, they constr ucted a photoion- 
izatio n model, using the Harrington code (e.g. iHarrington et alj 
Il982h . assuming a biconical geometry for a nebula seen almost 
pole-on. A two-dimensional representation of their model is shown 
in Figure 10 of C87. They had first tried a spherical model, deriving 



the radial hydrogen density distribution from the average radial in- 
tensity profile in the H/3 map they used, which had been obtained at 
the Boyden Observatory, South Africa, in 1973 by Drs. K. Reay and 
S. P. Worswick. The Zanstra temperature they derived for the cen- 
tral star was 117, 000 K, corresponding to a luminosity of 4900 1/0, 
for an adopted distance of 1.5 kpc. This model did not succeed 
in reproducing the line strengths from some of the high ionization 
species observed, such as, for example, Ne V, O V and O IV. They 
interpreted this as an indication that the nebula could be optically 
thin in some directions as seen from the star; this led to the for- 
mulation of a biconical model. The presence of an optically thin 
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Table 1. Input parameters for the NGC 3918 biconical model. 



position 



wavelength 



Figure 1. Model image (left) and lo ng axis sy n thetic slit spectrum (right) 
from the hydrodynamical models of ' Mellema il99d) . Figure adapted by 
C99 from Mellema's original thesis work. 
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phase required an upward correction to the original Zanstra tem- 
perature and, in the final model, they adopted the ionizing spec- 
trum described by a non-LTE model atmosphere for a central star 
having an effective temperature of 140,000 K, a surface gravity of 
log g = 6.5 and and a photospheric He/H = 0.10 by number. 
T his model atmosphere w as calculated by C87 using the program 
of lMihalas & AuedJ 19751) . The resulting nebular model seemed to 
reproduce the main spectroscopic features observed. 

Present observational data, such as for example images of 
NGC 3918 taken by the HST (see Figure [Sj and the echellograms 
obtained in several emission lines by iCorradi et al.. (J999, from 
now on C99) (Figure [Tot , show, however, that a biconical model 
is inconsistent with the spatio-kinematical structure of this plane- 
tary nebula. C99 presented an analysis of optical images and high 
resolution long slit spectra of three planetary nebulae, including 
NGC 3918. They concluded that the large scale structure of this 
object consists of a bright inner shell of roughly elliptical shape, 
from which two fainter protrusions extend in the polar directions, 
giving what was described in their paper as an overall spindle-like 
appearance. This inner shell, which has a size of 12" x 20", mea- 
sured along its major and minor axes, is surrounded by an ellip- 
tical shell with a diameter of 16". From the images and from the 
long slit spectra they obtained the basic kinematical, geometrical 
and orientation parameters of the inner shell. They adopted a hy- 
drodynamical model by G. Mellema (1996) to reproduce, at least 
qualitatively, the observations. The model which gave the best fit 
was one, from Mellema's set, which posits an initial mass distri- 
bution strongly enhanced in the equatorial region, with a density 
contrast between the equatorial and the polar regions as large as 
10 (see Figure0. The effects of an expanding shock driven by a 
strong wind would give a spindle-like structure similar to the one 
observed in the inner shell of NGC 3918. C99 derived the inclina- 
tion and the kinematical parameters of the inner shell by using the 
spatio-kinematical model of iSolf & Ulrich. l,.1985i) (which was also 
used to obtain Ha position-velocity diagrams, as well as the shape 
of the inner shell) to match the observational data. Their final model 
still showed some deviation from the observations, particularly in 
the long-axis position-velocity plot, which they attributed to sim- 
plified assumptions in the spatio-kinematical model. Photoioniza- 
tion calculations, however, were not carried out in their work, and 
therefore no comparison with the observed spectrum was available 
to them. 

Given the large amount of observational data available for 
this object and the existence of the two different models described 
above, NGC 3918 seemed an excellent candidate for a detailed 
three-dimensional photoionization study using the Mocassin code 



described by 'Ercolano| i2002h and 'Ercol ano et alj i20o3)- Three 
photoionization models were constructed, using different density 
distribution descriptions, in order to try to reproduce the main spec- 
troscopic features, as well as the projected maps published by C99. 
The models used a 23x23x23 Cartesian grid, with the ionizing 
source being placed in a comer in order to utilize the symmetry 
of the geometries used. The complete model nebulae were, there- 
fore, contained in 45x45x45 cubic grids. All the grid cells have 
the same size. Velocity fields were then applied to the final con- 
verged grids in order to produce position- velocity diagrams to com- 
pare with the observations. In Section 2 the biconical density dis- 
tribution model of C87 is described and the results obtained for it 
with Mocassin are presented and compared with those of C87, who 
used the one-dimensional Harrington code. A consistency test for 
the diffuse radiation field is also carried out in this section, in or- 
der to investigate the effects of discontinuities in the diffuse field 
transport in one-dimensional codes. The results obtained from the 
spindle-like models of NGC 3918 are presented in Section 3. Dis- 
crepancies were found between the predictions of all the models 
and the ISO SWS measurements of the infrared fine-structure lines; 
the possible reasons of these discrepancies are discussed in Sec- 
tion 4. 



2 THE BICONICAL DENSITY DISTRIBUTION MODEL 

The first nebular model to be run was the pole-on biconical distribu- 
tion described by C87. This consisted of an optically thick biconical 
structure embedded into an optically thin sphere. The bicones are 
aligned along the polar axis of the nebula, which describes an an- 
gle of 15° with the line of sight. Exactly the same parameters were 
used for the central star and, apart from helium, argon and iron, the 
same nebular abundances as in their model; these are summarised 
in Table □ The main reason for using exactly the same parame- 
ters is to be able to distinguish between the causes of any discrep- 
ancies between our model results and the original C87 model re- 
sults. However, in the case of helium the empirical value derived 
by C87 was used instead of their model value, since it was com- 
mented in their paper that this is more reliable. In the case of ar- 
gon and iron their empirical abundances also had to be used, since 
these elements were not included in their code (that of J.P. Harring- 
ton) and, therefore, photoionization model results for them did not 
exist. The biconical radial density law that was adopted by C87, 
and which was mapped onto the three-dimensional Cartesian grid 
used in the current model, was constructed by them by deconvolv- 
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Table 3. Diagnostic ratios for the electron density, N^, and electron temper- 
ature, Te (C87 biconical model). 
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2 


III 
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0.013 


0.013 


0.013 


1 


Ne IV 
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0.049 


0.049 


0.047 


2 


Ne V 


1575/3426 


0.0104 


0.012 


0.007: 


2 



Moc. = Mocassin results 

" In order of increasing ionization potential. 

^ (7321-i-7332)/(3726+3729) 

" (4068+4076)/ (67 1 7+673 1 ) 

references: l: lTsamis et al]<2003ft : 2: IClegg etai]<1987ft 



ing the observed nebular surface brightness, using the formalism of 
iHarrington & FeibelmanI ^S^), and has the form 



Nh{R) = iVo {R/RiT exp[-0.5(i?/i?i) 



(1) 



where Nh (R) is the hydrogen number density at distance R from 
the centre of the nebula, A'o is a normalization constant and Ri 
and n are parameters which may be varied. The values used for 
these parameters are also given in Table Q together with the 
half angle of the density enhanced cones, the density enhancement 
factor between the two regions, /, and the axial ratio, a/b. The 
optical depths at 1 . Ryd, calculated at the outer edge of the nebula in 
the polar (optically thick) and equatorial (optically thin) directions 
are 29.4 and 0.94, respectively. 



2.1 The Biconical Density Distribution Model: Results 

Following the C87 analysis, emission line fluxes are given relative 
to the total dereddened Hf3 flux integrated over the entire nebula, 
on a scale where H/S is equal to 100. The contributions from each 
sector and the total relative fluxes are listed in Table |2| together 
with C87's results for both sectors, and the total observed values. 
The last column of the table contains the references for the obser- 
vational data used. The fractional ionic abundances were also cal- 
culated for the thin and the thick regions and are given in Table|4] 
The mean temperatures weighted by the ionic abundances are cal- 
culated and given in Table|5| The definitions for both the fractional 
ionic abundances, f{Ni\), and for the mean temperatures, T{Nij, 
are given by Harrington et al. 1 1982) and are repeated here for con- 
venience: 



and 



/ NcNiidV 



(2) 



(3) 



where the average density of an ion over the nebular volume is 

r N^NijdV 

defined as < A'^i; > — ''—p ; A^c and Ni\ are the electron 

/ N^dV 

density and the the density of the ion j of element i, respectively, 
and the summation is over all k ionization stages; Tc is the electron 
temperature and dV is the volume element. 

The fractional ionic abundances and the mean temperatures 
weighted by ionic species were also presented by C87 (Tables 18 
and 19 in their paper), these are of great use in this work as they 
provide important parameters for the comparison of the Mocassin 
model to that of C8 7. The observed va lue of the integrated H/3 flux 
was obtained from 'Cahn et alj i 19921) and it was dereddened by 
adopting an interstellar reddening constant, c(H/?), equal to 0.40. 
This value, which is within the error bars of t he value estimated by 
C87, c{HI3) = 0.43 ± 0.05, was derived bv lTsamis et al. ( 20q|) 
and has been adopted in this work because it was derived from ob- 
servations of the optical Balmer emission lines, also published by 
[jsamis et al. ( 2003) and listed in our Table|2|(as reference a). It was 
not necessary to adjust the infrared line fluxes as these are not af- 
fected significantl y by interstellar red dening. The optical observa- 
tions presented bv lTsamis et alj i2003h were obtained by uniformly 
scanning a long-slit across the nebula, which method, when com- 
bined withthejntal H/3 flux measured using a large entrance aper- 
ture fe.g. lCahn et alJll992l) . yields the total fluxes from the whole 
nebula for all emission lines detected in the optical spectrum. These 
fluxes can thus be directly compared to those measured with space- 
borne facilities using large apertures, such as the lUE in the ultra- 
violet and ISO in the infrared. 

Most of the results shown in Table|2|are in reasonable agree- 
ment with the results published by C87, also reported in Table|2| 
For both models, however, the doubly ionized helium recombina- 
tion line at 4686 A is too strong, while the singly ionized helium 
triplet line at 5876 A is too weak, when compared with the obser- 
vations. This result is consistent with the ionic fractions listed in 
Table 18 of C87 and in Table |4] of this work for these ions. It is 
worth noting, however, that the values of the relative fluxes that 
C87 observed for these lines are 45.0 for He II4686 A and 9.5 for 
He I 5876A, on a scale where H/3 — 100. These values are closer 
to their model values of 48. and 7.8 than are the more recent obser- 
vations of 41.7 and 10.8 bv lTsamis et~ai]j2003h quoted in Tabled 
the former being the values they had tried to reproduce in their 
model. A better agreement with the more recent observations can 
be obtained by using a different model atmosphere for the ionizing 
source. As has already been mentioned, a non-LTE model atmo- 
sphere for a central star with an effective temperature of 140,000 K, 
log(? = 6.5 and a photospheric He/H = 0.10 by number was 
used by C87 and here; however by increasing log g, the depth of the 
doubly ionized helium jump is increased. This will result in fewer 
He^ ionizing photons and, therefore, lower He'^"'" and higher He^ 
nebular ionic fractions. The values obtained by using a non-LTE 
model for a central star of the same luminosity, effective temper- 
ature and He/H ratio, but with log g = 7.5, are 9.2 and 48.0 for 
the He I 5876A and He II 4686A line fluxes (relative to H/3, on a 
scale where H/3 = 100), respectively. These values are closer to the 
observations than those obtained with logg = 6.5 and reported in 
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Table 2. Observed and predicted emission lines fluxes for NGC 3918 and predictions for the C87 biconical models. 
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" Attenua ted by dust absorpti on. 

■ references: ai lTsamis et ^l2003l) : b: ICahn et al]<1992h : c: lLiu et aljhoOll) : d: lBoweill200ll) : e: C87; f: |Pottasch et'ai]<198d) . 

Table|2| One final point to be made about the nebular helium results the difference between the results we obtained for the He II 4686 A 
is that, as was mentioned at the start of this section, the nebular he- and He I 5876 A line fluxes and those presented by C87. 
lium abundance used in this work is the empirical value calculated 
by C87, He/H=0. 107 by number, which is slightly higher than the 

value they used in their model, He/H=0.10; this also contributes to Large discrepancies also exist between the two models' pre- 

dictions and the ISO SWS fluxes measured by Bower (2001 - ref- 
erence d in Table |2}. This discrepancy is mainly due to a pointing 
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Table 4. Nebular averaged fractional ionic abundances for NGC 3918 (Mocassin: C87 biconical model) 
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0.039 


0.214 
0.348 


0.081 
0.438 


0.022 
0.173 


Si 


0.21(-3) 
0.14(-6) 


0.325 
0.62(-3) 


0.195 
0.015 


0.256 
0.092 


0.127 
0.322 


0.078 
0.425 


0.019 
0.144 


Mg 


0.32(-2) 
0.61 (-6) 


0.112 
0.59(-4) 


0.532 
0.34(-2) 


0.119 
0.056 


0.118 
0.298 


0.098 
0.507 


0.018 
0.136 



For each element the upper row is for the optically thick phase and the lower row is for the optically thin phase. 



error which affected all ISO observations of NGC 3918; this is dis- 
cussed in more detail in Section 4 and Appendix C. 

With regard to the other emission line fluxes listed in Ta- 
ble |2| some discrepancies between the two models do exist, but 
they are, in the majority of cases within 30% and the overall ion- 
ization structure of the nebula coincides. In particular, large dis- 
crepancies were found between the current work and that of C87 in 
the predictions of the [Ne V] fine- structure line fluxes, which will 
be discussed in Appendix B. Considering the significant progress in 
atomic physics made in the past fifteen years, on the one hand, and 
the completely different treatment of the radiative transfer adopted 
by the two codes, on the other, that there are some differences in the 
results obtained is hardly surprising. Neither the Opacity Project 
^Verner et al. 1996) nor the Iron Project (Hummer et al. 1993) ex- 
isted in 1987 and the rates for other important processes have 
changed since then. 

Mocassin's predictions for electron density and temperature 
diagnostic line ratios can be compared to those obtained by C87, 
and to the observed ratios, in Table|3| There is very good agreement 
between the two models' predictions of line ratios, which is to be 
expected, since the density distributions used were the same. The 
small differences observed are reflected in the predicted emission 
line fluxes, listed in Table|2| and, as has already been discussed, can 
probably be mostly assigned to differences in the atomic data sets. 
The diffuse radiation field, however, might also have some effect 
on the ionization structure, especially at the edges between the two 
sectors and in the outer regions of the nebula. 



2.1.1 Resonance lines 

There are some large discrepancies between the predictions of both 
the Mocassin and C87 models and the observations for the fluxes 
of the C IV 1550 A and Mg II 2800 A resonance lines. The reason 
for this is that these lines are attenuated by nebular dust absorption, 
which was not accounted for in this work, or by C87. Moreover, 
Mg II 2800 A also suffers from absorption by interstellar Mg"'" and 
this is probably the main reason for the almost total absence of this 
line in the observed spectrum. Resonance lines must travel a much 
longer effective distance to escape from a nebula than do forbid- 
den lines or the continuum radiation, because of the high optical 
depth in them; this translates into an enhanced probability of ab- 
sorption by dust. The N V 1240 A resonance doublet is also affected 
by this process, although in a minor way. The Si IV resonance line 
doublet at 1393 A and 1401 A can also be subject to attenuation 
by dust absorption, however this line appears to be only slightly 
overestimated by the Mocassin model and is in agreement with the 
observations in C87. Harrington, Monk & Clegg published a dust 
model for NGC 3918 which confirms that the inten sities of the reso- 
nance lines are indeed reduced. In their dust model lHarrington et alj 
(1988) followed the same method as in the C87 photoionization 
model for this object; treating the thin and the thick sectors as two 
separate, spherically symmetric objects and then combining the re- 
sults, weighted by the respective volume of each sector, i.e. 29% 
for the thick sector and 7 1 % for the thin sector. In their model they 
found that the total flux of the C IV 1550 A doublet, relative to H/3 = 
100, was reduced from 1435 to 950, with 45% of photons destroyed 
in the optically thick sector and 21% in the optically thin sector. 
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Table 5. Mean temperatures (K) weighted by ionic species for NGC 3918 (Mocassin: C87 biconical model). 



Ion 

Element I n m TV V VI vff 



H 


11,386 
15,881 


12,925 
16,326 












He 


11,247 
15,495 


12,044 
15,722 


14,703 
16,344 










C 


11,582 
15,300 


11,585 
15,459 


12,116 
15,603 


13,133 
15,861 


15,858 
16,526 






N 


10,339 
15,348 


11,522 
15,506 


12,182 
15,631 


13,186 
15,842 


15,435 
16,234 


17,312 
17.217 







10,834 
14,874 


11,618 
15,298 


12,221 
15,509 


14,507 
15,887 


15,680 
16,260 


17,081 
17,020 


20,296 
19,477 


Ne 


10,136 
14,990 


11,636 
15,243 


12,156 
15,467 


14,340 
15,859 


15,561 
16,269 


17,360 
17,301 


22,911 
21,769 


S 


11,549 
15,016 


11,412 
15,115 


11,845 
15,353 


12,409 
15,603 


13,887 
15,910 


15,5584 
16,285 


17,637 
17,437 


Si 


11,817 
15,362 


11,794 
15,398 


12,349 
15,561 


12,875 
15,677 


14,333 
15,949 


15,543 
16,364 


17,684 
17,563 


Mg 


11,747 
14,517 


11,690 
15,530 


12,013 
14,963 


13,491 
15,632 


14,605 
15,919 


15,547 
16,329 


17,563 
17,526 



For each element the upper row is for the optically thick phase and the lower row is for the optically thin phase. 



Their attenuated flux was still 85% too strong compared to the ob- 
served value that they used (F(1550) = 512) and 107% too strong 
compared to the more rece ntly re-calibrated an d extracted lUE 
archive spectra presented by Tsamis et al. ('2003') (F(1550) = 458). 
[Harrington et al. 1 1988) did not give a final answer as to whether 
the cause of this discrepancy lies in the photoionization model pro- 
ducing too much C IV line radiation or in the dust model being 
too ineffective in destroying it, however they speculate that if the 
grain opacity near 1550 A was greater in their dust model, more 
line photons would be destroyed hence improving the fit. The flux 
of N V 1240 A, relative to H/3 = 100, is reduced from 66 to 52, 
with 25% of photons being destroyed in the thick sector and 20% 
in the thin sector. The value measured by Tsamis et al. (2003) for 
this line is 41.4. In the case of Si IV 1400 A dust actually makes 
the fit worse, by attenuating the flux, relative to H/3 = 100, from 8.2 
to 5.7, whereas the observed value is 9.0 (C87). This is, however, 
a weak line and the observed flux might not be accurate. Finally, 
[Harrington et al. ( 1988) predicted that the flux of Mg II 2800A is 
also reduced by dust attenuation, from 43 to 29. However, as men- 
tioned above, this line shows evidence of interstellar absorption, 
so it cannot be meas ured properly. Aside from this, the dusty pho- 
toionization model of lHarrington et all il988h is nearly identical to 
the one presented by CSV. 



2.1.2 The ionization structure 

The general ionization structures of the C87 biconical nebula cal- 
culated by the two codes, and shown in Table 18 of C87 and in 



Table|4]of this work, are in satisfactory agreement, with small dif- 
ferences which, again, can mainly be attributed to improvements 
in the atomic data and the different treatments of the radiative 
transfer. In fact, although both codes treat the diffuse field exactly, 
the C87 model could not allow for transfer of diffuse radiation 
between the thick and the thin sectors. As show n by the scatter 
of values obtained for the benchmark problems jPeauignot et alj 
[200 1; Ercolano 2002; Ercolano et al. 200^, even amongst one- 
dimensional classical codes it is usual for different codes to re- 
turn slightly different predictions for the same problem, even in the 
extremely simplified case of homogeneous, spherically symmetric 
nebulae. Figure|2|shows plots of the ionization structure of oxygen 
(middle panels) and carbon (bottom panels) as a function of radius 
along the optically thin polar direction (left panels) and along the 
optically thick equatorial direction (right panels). The ionization 
structure of helium is also plotted in each panel. The electron den- 
sity distribution is also plotted as a function of radius for the opti- 
cally thin polar direction (top left panel) and for the optically thick 
equatorial direction (top right panels). The effects of the charge 
exchange reaction 0'^"'"-i-H'^— >0^^-l-H"'" are evident, particularly in 
the optically thick region, where /(He^+) is about 10% larger than 
the sum of/(0^+)+/(0'*+)-t-/(0^+)-(-... 



2.1.3 The temperature structure 

The mean temperatures weighted by ionic species calculated by 
Mocassin for the C87 biconical model are listed in Table js] these 
can be compared to those calculated by C87, as listed in Table 19 of 



3-D photoionization modelling ofNGC 3918 7 




5,0X10'" i.Oxlo" I.5XI0'' 2.OKI0'' 2,5x10" 
radijs (cm) 





5.0x1d'^ 1.0x10''' 1.5x1g'' 

rodijs (cm^ 


.Oxio''' 2.5x1o''' 3.0x-o'^ 


r\ 

\ 06+ 


1 , , , , 1 , 


\ 

He2+ \ ^^^^ 


03+ 








02+ . 
==^^^^-~~^^^-=^^ ■ 


5,0x10'^ -.om^^ 1,5x10''' 2, 

rodijs (cm' 


xlO^^ 2,5x10'^ 3,0x10'^ 




Figure 2. Mocassin ionic fractions for oxygen (middle panels) and carbon (bottom panels) along the polar (left panels) and equatorial (right panels) directions 
for the C87 biconical model for NGC 3918. The ionization structure of helium is also plotted in each panel. The top panels show the electron density plotted 
as a function of radius along the polar (top left panel) and equatorial (top right panel) directions. 



their paper. The most obvious difference between the two sets of re- 
sults is that Mocassin's mean temperatures are systematically lower 
than Clegg et al.'s. This effect is even more pronounced for the op- 
tically thin region. First of all, the higher collision strengths for the 
[Ne v1 fine structure line s, as calculated bv lAggarwalHl983l) and by 
iLermon & Burkd il99lh . and used in this work (see Appendix B), 
will make [Ne V] an important coolant. Furthermore, the fact that 
in the optically thin region the fractional ionic abundance of Ne*+ 
is higher could also explain why this region seems to be more af- 



fected by the cooling. It is also worth noticing that the mean ionic 
temperatures predicted by C87 show very little variation for the 
neutral species in each sector. These temperatures are all concen- 
trated around 12, 065 ± 35 K for the neutral species in the optically 
thick sector, and around 16, 120 ± 210 K for the optically thin sec- 
tor (if re(ff°) = 16, 330 K is excluded then the rest of the neutral 
ion temperatures in the optically thin sector are grouped around the 
value 16,005±95 K). The scatter gets larger for the higher ioniza- 
tion species. Mocassin's results, on the other hand, show a larger 
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scatter, even for the mean temperatures of the neutral species, with 
their temperatures grouped around 10960 ± 896 K for the opti- 
cally thick region and 15105 ± 682 K for the optically thin region. 
One probable explanation for this effect is the fact that, in C87's 
work, the biconical geometry is reproduced by combining two sep- 
arate spherically symmetric models, hence taking no account of the 
interaction between the diffuse fields from the two phases, which 
could lead to fluctuations in local temperatures, particularly at the 
boundaries between the two regions. This is discussed in more de- 
tail in the next subsection. 



2.2 Diffuse Radiation Field Consistency Test 

C87 combined two appropriately weighted, spherically symmet- 
ric models to reproduce the biconical geometry. A consequence of 
this approach was that the diffuse radiation field was not treated 
self-consistently, especially near the boundaries of the two zones. 
The effects of the diffuse radiation field are expected to be small, 
but they could account for some of the differences between Moc- 
cassin's results and those obtained with the spherically symmet- 
ric Harrington code. In order to investigate the magnitude of these 
effects, the biconical distribution was modelled again using Moc- 
cassin, but this time each sector was treated as if it were spherically 
symmetric. The results from each sector were then weighted by the 
volume of the same sector and added together to yield the total flux 
relative to the total H/3 flux from the whole nebula. From now on 
this model will be referred to as the composite model. 

Table |6| lists the fluxes of the emission lines relative to H/3, 
on a scale where H/3 is equal to 100, for each sector For ease of 
comparison, the results from the previous, self-consistent model 
are also given in the same table. Tables and |8| list the fractional 
ionic abundances and the weighted mean electron temperatures for 
individual ionic species for the composite model. It becomes ap- 
parent from inspection of Table |6| that all of the collisional lines 
coming from lower ionization species are much weaker in the op- 
tically thick region of the composite model than in the same region 
of the self-consistent model. This is consistent with the fractional 
ionic abundances of low ionization species being lower in the opti- 
cally thick region of the composite model (Table[8}, hence yielding 
lower line fluxes. The helium recombination line intensities cal- 
culated for the optically thick region of the composite model also 
indicate a lower abundance of He" and He^^, but a higher He''" 
abundance, which is confirmed by the results in Table0 The opti- 
cally thin region results do not seem to be affected as much and they 
are almost identical for the composite and for the self-consistent 
model. There are two main causes for the discrepancies observed 
between the results for the optically thick regions of the composite 
and self-consistent models. First of all, the diffuse radiation field 
coming from the optically thin region of the self-consistent model 
is weaker than the diffuse field coming from the same region in 
the optically thick spherically symmetric nebula of the composite- 
model, since the diffuse field is approximately proportional to the 
square of the density and the densities in the two regions of the 
spherically symmetric optically thick composite model are, obvi- 
ously, the same and are higher than the density in the optically thin 
region of the self-consistent model. This means that more diffuse 
photons are available to ionize the low ionization potential species 
in the composite model than there are in the self-consistent one. 
Secondly, the shape of the diffuse continuum at the edges between 
the two regions is different in the two models. In the self-consistent 
model the diffuse photons coming from the optically thin region 



Table 6. Predicted emission lines fluxes for NGC 3918 (composite versus 
self-consistent model.) 



Line Composite Self-Consistent 

Thin Thick Tot. Thin Thick Tot. 



will have higher frequencies, on average, than those coming from 
the same region in the optically thick spherically symmetric nebula 
of the composite model; this is due to the fact that in the optically 
thin medium the gas is in a iiigher ionization state. The photoion- 
ization cross-sections of the lower ionization species in the opti- 
cally thick region of the self-consistent model will be, in general, 
quite small at the higher frequencies typical of the diffuse photons 



cm^ s 

H/3 4861 
He I 5876 
He II 4686 
C II] 2326 
C II 4267 
[C II] 157.74/im 
C III] 1908 
C III 4648 
C III] 2297 
Civ 1550 
C IV 4658 
[Nl] 5199 
[N II] 6584 
Nlll] 1751 
[N III] 57.3 iim 
N iv] 1486 
N V 1240 
[O l] 6300 
[Ol] 63.12 /^m 
[O l] 145.5 /xm 
[O II] 3727 
[Olll] 51.8 ^tm 
[O III] 88.4 ^tm 
[0 III] 5007 
[0 III] 4363 
[0 IV] 25.9 Aim 
O iv] 1402 
O IV 1342 
[Ne II] 12.8 Atm 
[Ne III] 15.6 /xm 
[Ne III] 3869 
[Ne IV] 2423 
[Ne v] 3426 
[Ne v] 14.3 Aim 
[Ne v] 24.3 Aim 
[Ne Vl] 7.6 Aim 
Mg I] 4565 
Mg II 2800 
[Mg iv] 4.49 Atm 
[Mg v] 2783 
[Mg v] 5.6 Aim 
Si III] 1887 
Si IV 1399 
[S II] 4070 
[S III] 18.7 Aim 
[S III] 33.6 Aim 
[S III] 6312 
[S iv] 10.5Aim 



0.45 


1.90 


100 


100 


0.07 


8.7 


23 


28 


0.22 


34 


0.004 


0.31 


.0005 


.12 


33 


551 


0.10 


0.56 


1.05 


6.6 


607 


968 


0.24 


0.16 


0.0 


0.006 


0.06 


27 


2.0 


22 


0.60 


14 


19 


26 


44 


21 


0.0 


0.18 


0.0 


0.011 


0.0 


.0008 


0.08 


47.5 


2.4 


101 


0.75 


19. 


35 


1624 


0.78 


20 


186 


157 


16 


14 


1.3 


0.84 


0.001 


0.62 


1.6 


130 


3.7 


161 


36 


41 


76 


44 


176 


110 


154 


65 


49 


21 


0.0 


0.34 


0.20 


57 


0.39 


3.1 


2.2 


2.4 


4.4 


6.6 


1.1 


14 


2.9 


7.6 


0. 


1.3 


0.03 


32 


0.020 


14 


0.006 


3.2 


6.8 


223 



2.35 


0.45 


100 


100 


8.8 


0.08 


51 


24 


34 


0.21 


0.31 


0.004 


.12 


0.0 


584 


33 


0.66 


0.10 


7.7 


1.11 


1575 


596 


0.40 


0.25 


0.006 


0.0 


28 


0.06 


22 


2.0 


15 


0.64 


45 


19 


65 


42 


0.18 


0.0 


0.012 


0.0 


.0008 


0.0 


48 


0.08 


103 


2.8 


20 


0.91 


1659 


40 


21 


0.85 


143 


190 


29 


16 


2.1 


1.4 


0.64 


0.001 


133 


1.8 


165 


4.1 


77 


35 


120 


75 


286 


179 


219 


156 


70 


49 


0.38 


0.0 


57 


0.26 


3.5 


0.42 


4.6 


2.2 


11 


4.5 


15 


0.97 


11 


2.7 


1.3 


0.0 


32 


0.03 


14 


0.02 


3.2 


0.006 


235 


7.23 



1.88 


2.33 


100 


100 


8.5 


8.5 


30 


54 


57 


57 


0.32 


0.32 


0.35 


0.35 


577 


610 


0.50 


0.61 


5.9 


7.0 


964 


1560 


0.16 


0.40 


0.21 


0.21 


73 


73 


22 


24 


12 


13 


26 


45 


22 


64 


2.7 


2.7 


0.23 


0.23 


0.018 


0.018 


104 


104 


90 


93 


16 


17 


1543 


1583 


20 


20 


161 


352 


14 


30 


0.83 


2.2 


0.94 


0.94 


130 


132 


161 


165 


42 


77 


45 


119 


107 


286 


64 


220 


22 


70 


0.73 


0.73 


78 


78 


3.1 


3.5 


2.4 


4.6 


6.3 


11 


15 


16 


7.9 


11 


2.1 


2.1 


38 


38 


16 


16 


3.7 


3.7 


207 


214 
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Table 7. Fractional ionic abundances for NGC 3918 (Mocassin: composite biconical model) 



Ion 

Element I n m TV V VI vfT 



H 


0.61(-2) 
0.54(-3) 


0.994 
0.999 












He 


0.19(-2) 


0.697 
0.025 


0.301 
0.975 










C 


0.61(-4) 
0.56(-7) 


0.029 
0.20(-3) 


0.481 
0.027 


0.366 
0.254 


0.124 
0.718 






N 


0.97(-4) 
0.26(-7) 


0.025 
0.12(-3) 


0.476 
0.037 


0.399 
0.334 


0.068 
0.369 


0.032 
0.259 




O 


0.49(-3) 
0.16(-7) 


0.031 
0.72(-4) 


0.730 
0.034 


0.149 
0.360 


0.066 
0.403 


0.020 
0.167 


0.42(-2) 
0.035 


Ne 


0.38(-5) 
0.24(-8) 


0.76(-2) 
0.45(-4) 


0.747 
0.033 


0.126 
0.257 


0.098 
0.539 


0.020 
0.165 


0.61(-3) 
0.49(-2) 


S 


0.32(-5) 
0.16(-7) 


0.036 
0.14(-3) 


0.368 
0.015 


0.387 
0.126 


0.127 
0.275 


0.061 
0.413 


0.019 
0.170 


Si 


0.13(-3) 
0.19(-6) 


0.268 
0.66(-3) 


0.207 
0.016 


0.282 
0.094 


0.143 
0.319 


0.080 
0.425 


0.019 
0.144 


Mg 


0.17(-2) 
0.40(-5) 


0.082 
0.34(-3) 


0.638 
0.024 


0.133 
0.162 


0.066 
0.248 


0.065 
0.441 


0.014 
0.123 



For each element the upper row is for the optically thick phase and the lower row is for the optically thin phase. 



coming from the optically thin region and so these species will be 
less efficiently ionized by them. In the composite model, however, 
since the density distribution is spherically symmetric, the ioniza- 
tion structure is also the same in all radial directions; this means 
that the diffuse photons coming from adjacent regions will have fre- 
quencies close to the ionization thresholds of the lower ionization 
species which will, therefore, be ionized to a greater degree than 
in the self-consistent model. On the other hand, the optically thin 
region of the self-consistent model will not be affected by the lower 
energy diffuse photons coming from the optically thick region since 
the gas in the optically thin region is already in a higher ionization 
stage. The combination of the two effects discussed above explains 
why the abundances of lower ionization species are lower in the op- 
tically thick region of the composite model than in the same region 
of the self-consistent one. 

Another striking feature which emerges from running the 
composite spherically symmetric model is in the prediction of mean 
temperatures weighted by ionic species shown in Table|8| As men- 
tioned in the previous subsection, the mean ionic temperatures pre- 
dicted by C87 show very little variation amongst the neutral species 
in each sector, whereas the scatter of the results from Mocassin's 
self-consistent model is much larger. It was anticipated that a pos- 
sible explanation for this effect could be the fact that in Clegg et 
al.'s work the diffuse radiation field transfer was not self-consistent, 
particularly near the boundaries between the two density regions 
and that the interaction between the two phases' diffuse radiation 
fields could lead to the fluctuations in the kinetic temperatures ob- 
tained in Mocassin's self-consistent model. The results obtained 



from Mocassin's composite model indeed confirm this hypothesis 
- within each sector the fluctuations in the mean ionic tempera- 
tures weighted by ionic abundances for neutral species disappear 
and their distribution resembles that obtained by C87 although the 
Mocassin temperatures are still slightly lower, particularly in the 
thin sector. 

In conclusion, it is clear that the effects of the diffuse radia- 
tion field near the boundaries between the two phases are not neg- 
ligible, even for the simple biconical distribution used in this case. 
These effects are particularly important for the abundances of the 
lower ionization species. It is not clear, however, why the ionization 
structure calculated by C87 (their Table 18), who used a compos- 
ite model, bears more resemblance to Mocassin's self-consistent 
model results (see Table |4j, than to Mocassin's composite model 
results (see Table0. Unfortunately, due to the completely different 
approach to the radiative transfer used by the two codes and the 
fact that fifteen years have passed since C87's work, over which 
period the atomic database has evolved considerably, any further 
comparison between the two codes becomes difficult. 



3 A SPINDLE-LIKE DENSITY DISTRIBUTION MODEL 

FOR NGC 3918 

FigureQ from C99, shows a qualitative representation of a density 
distrib ution obtained from a hydrodynamical model by Mellema 
I* 199(1) that was adopted by C99 in their observational study of the 
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Table 8. Mean temperatures (K) weighted by ionic species for NGC 3918 (Mocassin, composite biconical model). 



Ion 



Element 



H 


12,060 
16,061 


12,880 
16,486 












He 


11,965 
15,844 


12,052 
16,001 


14,770 
16,498 










C 


11,964 
15,702 


11,985 
15,757 


12,108 
15,842 


12,931 
16,039 


15,661 
16,667 






N 


12,022 
15,704 


12,002 
15,761 


12,147 
15,842 


12,992 
16,006 


15,279 
16,369 


17,065 
17,360 







12,032 
15,678 


11,989 
15,749 


12,155 
15,838 


14,481 
16,036 


15,507 
16,388 


16,862 
17,153 


20,001 
19,615 


Ne 


11,958 
15,686 


11,990 
15,736 


12,146 
15,817 


14,301 
16,010 


15,394 
16,406 


17,122 
17,449 


22,502 
21,943 


S 


11,986 
15,675 


11,982 
15,689 


12,946 
15,742 


12,233 
15,850 


13,638 
16,075 


15,402 
16,419 


17,367 
17,573 


Si 


11,919 
15,674 


11,947 
15,704 


12,169 
15,773 


12,624 
15,858 


14,066 
16,100 


15,293 
16,518 


17,397 
17,744 


Mg 


11,919 
15,650 


11,929 
15,671 


11,957 
15,717 


13,408 
15,841 


14,490 
16,066 


15,372 
16,473 


17,322 
17,680 



For each element the upper row is for the optically thick phase and the lower row is for the optically thin phase. 



Table 9. Input shape and thickness parameters for the analytical fit to the 
density distribution of spindle-like model B (see Appendix A). 



ni [cm 
di [dcg] 



' eg 
' po 



18000 

25 

0.3 

0.2 
0.25 



A e [deg] 
Ar'" 

eg 
A ^out 
eg 

Ar™ 

A ^out 
po 



10 
0.06 
0.2 
0.0032 
0.5 



normalized to the outer radius of the nebula 



Table 10. Input parameters for the spindle-like models of NGC 3918. 





6900Lq 


N/H 


1.5x10" 


-4 


L*(B) 


578OL0 


O/H 


5.0x10" 


-4 


Rin 


0.0 pc 


Ne/H 


1.2x10" 


-4 


Rout 


0.106 pc 


Mg/H 


1.4x10" 


-4 




140,000 K 


Si/H 


1.0x10" 


-4 


logg 


7.5 


S/H 


1.6x10" 


-4 


He/H 


0.107 


Ar/H 


2.0x10" 


-4 


C/H 


8.0x10"* 


Fe/H 


3.7x10" 


-4 



The abundances are given by number, relative to H 




Figure 3. Three-dimensional iso-density plot of the dense torus component 
of spindle-like model A. 



spatio-kinematical properties of NGC 3918. The model chosen was 
one from a set that Dr. G. Mellema computed for his original thesis 
work (parameters B — 6 and A — 0.9 in his calculations, page 
126). 

The exact density distribution for the model used in the C99 
paper was kindly provided by G. Mellema, and a first attempt to 
construct a photoionization model for this distribution was carried 
out. The density distribution, however, had to be scaled up by a 
factor of 1.6 in order to match the observed H/3 flux. However, 
the integrated emission line spectrum still did not match the ob- 
servations, as the emission line strengths from the lower ionization 
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Table 13. Fractional ionic abundances for NGC 3918 (Spindle-like models of NGC 3918) 



Ion 



Element I 


n 


III 


IV 


V 


VI 


VII 


H 0.14(-1) 
0.85(-2) 


0.986 
0.991 












He 0.64(-2) 
0.32(-2) 


0.509 
0.513 


0.485 
0.483 










C 0.73(-4) 
0.96(-4) 


0.33(-l) 
0.23(-l) 


0.284 
0.258 


0.462 
0.509 


0.221 
0.210 






N 0.53(-2) 
0.23(-2) 


0.32(-l) 
0.25(-l) 


0.279 
0.256 


0.501 
0.544 


0.129 
0.125 


0.053 
0.047 




O 0.13(-1) 
0.71(-2) 


0.28(-l) 
0.24(-l) 


0.544 
0.562 


0.250 
0.257 


0.120 
0.112 


0.038 
0.033 


0.62(-2) 
0.52(-2) 


Ne 0.68(-3) 
0.12(-3) 


0.14(-1) 
0.66(-2) 


0.553 
0.560 


0.195 
0.207 


0.192 
0.186 


0.045 
0.039 


0.99(-3) 
0.92(-3) 


S 0.11 (-4) 
0.15(-4) 


0.21(-1) 
0.13(-1) 


0.116 
0.094 


0.295 
0.276 


0.386 
0.445 


0.149 
0.144 


0.032 
0.028 


Si 0.81 (-4) 
0.84(-4) 


0.133 
0.101 


0.128 
0.110 


0.247 
0.243 


0.211 
0.249 


0.211 
0.230 


0.070 
0.068 


Mg 0.15(-2) 
0.15(-2) 


0.59(-l) 
0.44(-l) 


0.351 
0.313 


0.123 
0.149 


0.179 
0.203 


0.225 
0.233 


0.060 
0.056 



For each element the upper row is for model A and the lower row is for model B 

behind it, and the latitudinal extension of the density enhancement 
was determined such that the model would return the observed H/3 
flux. This trade-off between the density enhancement factor and its 
spatial extension resulted in the creation of a thick equatorial waist. 
The three-dimensional iso-density plot (obtained for a high density 
threshold) of the final density distribution used for this model is 
shown in Figure|3| This model will be referred to as model A. 

A second attempt to reconstruct the density distribution of 
this object was to seek an analytical expression which could re- 
produce the structure shown in Figure Q and in the HST and NTT 
images published by C99. The parameters of the central star and 
the adopted nebular abundances were the same as in the bicon- 
ical model. The density distribution was obtained by combining 
a sphere with an ellipsoid, using a set of analytical epressions as 
described in Appendix A. The best fit to the observed integrated 
spectrum and to the images of NGC 3918 was obtained using the 
parameters listed in Table|9| Figure|4]shows the three-dimensional 
iso-density plot corresponding to this combination. In this work the 
spindle-like model constructed by using the analytical description 
will be referred to as model B. 

C99's study primarily focused on the kinematics of NGC 3918 
and therefore they did not propose any new values for the cen- 
tral star parameters or for the nebular elemental abundances. Al- 
though derived for a different density distribution, the parame- 
ters derived by C87 were used by us for these models too. As 
given by C87, the dereddened stellar flux at 1300 A is equal to 
2.4x 10~^^ergscm~^s^^A^^. For an assumed distance of 1.5 kpc 




Figure 4. Three-dimensional iso-density plot of the spindle-like density dis- 
tribution for model B.The density cut was chosen in order to show the shape 
of the inner shell of the nebula. 



species were heavily underestimated. From the images, however, 
it is evident that NGC 3918 has an inner shell of denser material 
which screens the gas behind it from direct starlight, thus allow- 
ing the existence of a lower ionization region. In Mellema's sim- 
ulations such shells can only be made by considering an evolving 
central star and having a slow H ionization front. This, however, 
was not the case for the model used in C99 analysis. In order to 
simulate such a dense inner shell, the density was artificially raised 
in the equatorial region to form a thick torus between the two lobes. 
The torus had to be thick enough to effectively screen the material 
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Table 14. Mean temperatures (K) weighted by ionic species for NGC 3918 (Spindle-like models of NGC 3918) 



Ion 

Element I n m TV V VI vn" 



H 


10,142 
10,763 


13,739 
13,593 












He 


9,795 
10,519 


12,430 
12,457 


15,061 
14,770 










C 


10,683 
10,945 


10,715 
1 1,193 


12,521 
12,484 


13,733 
13,518 


15,539 
15,278 






N 


9,254 
9,805 


10,853 
11,239 


12,653 
12,606 


13.805 
13,594 


15,381 
15,085 


16,088 
15,924 







9,606 
10,188 


11,256 
11,518 


12,729 
12,689 


14,939 
14,656 


15,524 
15,288 


16,020 
15,792 


17,278 
17,839 


Ne 


9,093 
9,663 


10,202 
11,186 


12,603 
12,574 


14,771 
14,464 


15,422 
15,149 


16,004 
15,785 


17,881 
19,909 


S 


10,487 
10,872 


10,207 
10,796 


11,871 
11,840 


12,801 
12,712 


14,221 
13,875 


15,434 
15,133 


16,218 
16,117 


Si 


11,636 
11,628 


11,557 
11,652 


12,579 
12,443 


13,117 
12,968 


14,308 
13,935 


15,081 
14,655 


15,736 
15,381 


Mg 


11,403 
11,446 


11,189 
11,456 


12,197 
12,106 


13,730 
13,338 


14,759 
14,384 


15,250 
14,902 


15,812 
15,529 



For each element the upper row is for model A and the lower row is for model B 



and a stellar effective temperature of Te//= 150,000 K, this flux cor- 
responds to a stellar luminosity, L,, of 5780 L0. 

The stellar luminosity in the C87 model, how- 
ever, was finally adjusted by them to a slightly higher 
value of F.(1300A) = 2.9xlO"^^ergscm"^s"^A"\ cor- 
responding to L, = 6900 Lq. In this work both values 
for F,(1300A) were used, in turn, for both models A 
and B. Finally, F*(1300A) = 2.9x lO^^^ergscm^^s^^A"^ 
(L« = 6900 Lq) was chosen for model A and 
F.(1300A) = 2.4xlO~^^ergscm-2s-^A-i (L* = 578OL0) 
was chosen for model B, as these produced the best results. 
Table [Tol summarises the nebular and stellar parameters used for 
the models described in this section. 



3.1 Spindle-like Density Distribution Models: Results 

Table ll ll lists the predicted emission line fluxes integrated over the 
nebular volume for the spindle-like models A and B for NGC 3918; 
for comparison, the observed values are also reported in the last 
column of this table. The observed spectrum is reasonably repro- 
duced by both models. The H/3 fluxes predicted by models A and 
B are slightly different, but they are both within 10% of the ob- 
served value. The emission line predictions for helium seem to in- 
dicate that in both models the abundance of He^^ is slightly over- 
estimated, hence leaving too little He^. However both results are 
comparable to those obtained using the biconical model and indeed 
to the C87 results presented in Section lzTI The major discrepancies 
between the two models and the observations are the same as for 



the biconical model case. The C IV 1549 A and the Mg II 2800 A 
resonance doublets are overestimated by large factors by the Mo- 
cassin models, which do not take dust absorption into account. The 
Si IV 1393 A and 1401 A resonance line doublet, which should 
also suffer from absorption by dust, does not seem to be aifected . 
In fact if the 30% attenuation predicted by Harrineton et alj il988l) 
were to be applied to the flux of this doublet, the fit to the obser- 
vations would be made worse. However, as has already been dis- 
cussed in section |2. 1.11 Si IV 1393 A and 1401 A are quite weak 
and the flux measured might not be accurate. As before, large dis- 
crepancies with the observations, of factors of about ten, are also 
observed for the [Ne V] fine structure lines, namely [Ne V] 14.3 /im 
and [Ne V] 24.3/im; this is discussed in Appendices B and C. The 
large discrepancies found between the spindle-like model predic- 
tions for the infr ared line fluxes on the one hand, and the ISO SWS 
measurements of lBoweil i200ll) reported in Table ll ll on the other 
hand, are most probably due to the ISO SWS aperture having been 
offset from the centre of the nebula. This is discussed in more de- 
tail in Section|4]and Appendix C. In general, the integrated spectra 
obtained from the spindle-like models are in good agreement with 
those obtained from Mocassin's biconical model, and, therefore, 
comments made during the discussion of the latter (Section |2.H 
largely apply to the models A and B discussed in this section. The 
diagnostic line ratios calculated for models A and B are given in 
Table[T2l They are in satisfactory agreement with the observed ra- 
tios, also reported in the same table, with model B providing the 
slightly better fit. 

Tables fTsl and 1 1 41 list, respectively, the fractional ionic abun- 



3-D photoionization modelling ofNGC 3918 13 




5,0x10 ^ 1.0x10^^ 1.5x10^^ 2,Oxlo'^ 2.5x10^^ 5,0x10 ^ I.Oxlo'^ 1,5x10^^ 2,Ox1o'' 2,5x10^^ 

rodijs (cm) mdius (cm) 




rodijE (cm) rodijs (cm) 



Figure 5. Ionic fractions for oxygen (middle panels) and carbon (bottom panels) along the polar (left panels) and equatorial (right panels) directions. The 
ionization structure of helium is also plotted in each panel. The top panels show the electron density plotted as a function of radius along the polar (top left 
panel) and equatorial (top right panel) directions. Spindle-like model A for NGC 3918. 



dances, /(X'), and the mean electron temperatures weighted by 
ionic species, for the spindle-like distribution models of NGC 3918. 
For each element listed in the two tables, the upper entry is for 
model A and the lower entry is for model B. It is difficult to make a 
comparison with Mocassin's biconical model results or with those 
of C87, since both of the spindle-like models consist of a contin- 
uously varying gas density distribution where clear edges do not 
exist. The values of the ionic fractions reported in Table lOl were ob- 
tained from Equation J2}, and represent an average over the whole 



ionized volume of the local ionic fractions weighted by the local 
density. Figures|5|and|6|show the ionic fractions of oxygen (middle 
panels) and carbon (bottom panels) as a function of radius along 
the polar (left panels) and equatorial (right panels) directions for 
the spindle-like models A and B, respectively. The electron density 
distribution along the two directions is also shown in the top panels 
of Figures |5| and |6| for the polar (left panels) and equatorial (right 
panics) directions. It is clear from the figures that the ion distri- 
butions are very different along the two directions of each model, 
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O.B 
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G 
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0.2 
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Figure 6. Ionic fractions for oxygen (middle panels) and carbon (bottom panels) along the polar (left panels) and equatorial (right panels) directions. The 
ionization structure of helium is also plotted in each panel. The top panels show the electron density plotted as a function of radius along the polar (top left 
panel) and equatorial (top right panel) directions. Spindle-like model B for NGC 3918. 



according to the gas density structure. Integrated ionic fractions, 
/(X*), such as those hsted in Table llSl are very useful when dis- 
cussing the large-scale ionization structure of models with density 
variations, however plots of the type shown in Figures|5]and|6|can 
provide a better insight into the actual ionization topography of a 
model. 

As shown in Table fT4l the weighted mean temperatures cal- 
culated by Mocassin for model B (lower entries) are consistently 
lower than those calculated for model A (upper entries), apart from 



the third and fourth column values (for doubly and triply ionized 
species), which are, on the contrary, higher for model B. It is also 
evident that the scatter of values for the same column of ionic stage 
is larger for model B even for the lower ionization stages, while the 
scatter is very small for model A. The effects of density distribu- 
tions on the final values obtained for the mean weighted tempera- 
tures have already been discussed in Section l2^ It is clear that in 
this case the behaviour of the two models depends heavily on their 
density distributions. 
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Table 11. Observed and predicted emission lines fluxes for NGC 3918 
(spindle-like models). 



Table 12. Diagnostic ratios for the electron density, A^e , and electron tem- 
perature, Te (NGC 3918 spindle-like models). 





Predicted 






Line 


A 


B 


Observed 


ref'' 


H«/10""' 

^ cm'^ s 


2.19 


2.58 


2.34 


b 




100 


100 


100 




nc 1 J o / o 


8.15 


7.56 


10 81 






57.2 


53.6 


41.65 




P ttI 9'^9^^ 


34.6 


27.8 


31 24 




C II 4267 


0.23 


0.21 


0.50 




fr tt] 1 -sv 74/;rn 


0.22 


0.10 


0.21 




C tttI IQOS 

lllj 1 7UO 


543 


480 


492 5 




C III 4648 


0.89 


0.98 


0.42 




C tttI 97Q7 

111 / 


10.2 


11.3 


5- 




r TV 1 550° 


2267 


2220 


457.6 




r TV 4658 

1 V fUJO 


0.37 


0.37 


35 




fN t1 S 1 99 


1.03 


0.46 


0.63 




FN ttI fiSK4 


41.8 


33.9 


84 




K tttI 1 7S1 


22.9 


20.2 


zo. / 




[N III] 57.3 


7.41 


9.16 


13.5 




N Tv1 14S(i 

i'l 1 V 1 iT-OU 


61.4 


59.3 


46.4 




N V 1240^^ 


45.6 


38.7 


41 4 






3.95 


2.57 






\\J IJ UJ.IZ, ^111 


0.40 


0.26 


6 72 




rO t1 14S S //m 

IJ IH-J.J ^111 


.029 


.014 


23 




\Ci ttI "^797 

[w llj D 1 ^ 1 


42.7 


37.9 


94 5 




[0 III] 51.8 


79.4 


112 


88.6 




To tttI 9.9, 4 /;m 

W 111 00.*+ lAlll 


14.2 


22.2 


27 5 




\Ci tttI snn7 

V / 111 J UVJ / 


1760 


1789 


1657 




fO ml 4363 

\KJ All! 'T-JIJ-' 


25.0 


24.8 


21.6 




rn Tv1 9^ Q //m 


352 


365 


89 




iv] 1402 


34.9 


31.8 


48.9 




O IV 1342 


1.87 


1.72 


2- 




[Np TtI 19 S //m 

i'lC 11 1Z..O Lilll 


1.61 


0.75 


12" 


f 


TNTp TTt1 1 S ^^ //m 

i'lC 111 1-J.U tAlll 


123 


127 


46.3 


(J 


TNp ml 3869 

I'lC 111 -J\j\jy 


174 


173 


135.5 




\Ne ivl 2423 


89.2 


87.5 


132 




fNp vl '^49(1 

i'lC V .J'T'Z.U 


108 


99.1 


80 




[Ne v] 14.3 ^m 


267 


256 


19.6 


d 


[Ne v] 24.3 ^m 


170 


146 


23.1 


d 


[Ne Vl] 7.6 fim 


62.0 


54.4 






Mg I] 4565 


0.32 


0.33 


0.13 


a 


Mg II 2800"" 


47.2 


36.7 


^1 


e 


[Mg iv] 4.49 fj.m 


4.25 


5.10 






[Mg V] 2783 


4.84 


5.02 


5.6 


e 


[Mg v] 5.6 Aim 


12.5 


14.1 






Si III] 1887 


13.9 


11.2 


8 


e 


Si IV 1399'^ 


11.5 


10.2 


9 


e 


[S II] 4070 


1.90 


1.42 


2.69 


a 


[S III] 18.7 /im 


16.6 


14.3 


8.7 


d 


[S III] 33.6 fim 


5.21 


4.68 






[S III] 6312 


1.98 


1.55 


2.27 


a 


[S IV] 10.5Atm 


198 


213 


35.0 


d 



" Attenuated by dust absorption. 

^ refe ren ces: a: [Tsamis et alj l2003h: b: ICahn et alj <1992ft : 
1200 ll) : d: lBowejl200ll) : e: C87: f: IPottasch et alJ<1986ft . 



: iLiu et all 



Ion° 
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Ratio 


Ratio 


Ratio 


Ratio 


obs. 




(A) 


Moc. A 


Moc. B 


C87 


obs. 


ref.'' 




Mgl 


4562/4571 


0.38 


0.40 


0.42: 


0.28 


1 


Sll 


6716/6731 


0.65 


0.60 


0.70 


0.57 


1 


On 


3726/3729 


1.96 


2.02 


1.61 


1.98 


2 


On 


7325/3727'' 


0.11 


0.11 




0.0783 


1 


C III 


1906/1909 


1.25 


1.34 


1.36 


1.31 


2 


Si III 


1883/1892 


1.15 


1.31 


1.30 


1.37 


2 


Ne IV 


2421/2424 


1.09 


1.06 


0.93 


1.00 


2 


Niv 


1483/1486 


1.45 


1.47 


1.87 


1.85 


2 


IV 


1401/1404 


1.32 


1.36 


1.14 


1.50 


2 




Nil 


5755/6584 


0.021 


0.023 


0.025 


0.021 


1 


S II 


4073/6724^= 


0.200 


0.245 




0.339 


1 


C2+ 


4267/1908 


0.0004 


0.0004 


0.0005 


0.0008 


2 


S III 


6312/9532 


0.032 


0.032 


0.029 


0.032 


2 


III 


4363/5007 


0.014 


0.014 


0.013 


0.013 


1 


Ne IV 


1602/2423 


0.053 


0.050 


0.049 


0.047 


2 


Ne V 


1575/3426 


0.009 


0.009 


0.012 


0.007: 


2 



" In order of increasing ionization potential. 
(7321+7332)/(3726-i-3729) 

" (4068+4076)/ (67 17-1-673 1) 

<* references: L lXsamis et al]<2003ft : 2: C87 



4 THE INFRARED FINE-STRUCTURE LINE FLUX 

DISCREPANCIES 

The main reason for the discrepancy observed between the mod- 
els' predictions of the infrared line fluxes and the Bower (2001) 
measurements from the ISO SWS data for NGC 3918 is that the 
nebular coordinates used for the ISO observations were incorrect 
by approximately 14 arcsec (see Appendix C). A masking IDL 
routine ( Morisset et al. 2002), which makes use of the measured 
wavelength-dependent SWS beam profiles, was applied to the pro- 
jected maps obtained from the final converged model grids, in 
the six infrared fine-structure lines. This routine returns a mask 
which takes into account the wavelength-dependent aperture size 
and beam profile, and which can then be applied to a map of the ob- 
ject in the required emission line to predict the emergent flux that 
would be measured. The position of the aperture centre was off- 
set from the centre of the nebula by the amounts specified above, 
and the aperture was rotated by an angle of 80° East of North^ to 
correspond to the spacecraft roll-angle. Resultant predicted fluxes 
were obtained by convolution of the offset, rotated mask with the 
correctly oriented map of NGC 3918 in the required emission line. 



1 From the ISO archive, CINSTROLL = 350°, where CINSTROLL is the 
instrument's roll angle, defined as the angle anticlockwise between north 
and the spacecraft z-axis. Since the long axis of the SWS apertures is in 
the y-direction and the short axis is in the z-direction, then P.A. = CIN- 
STROLL -H 90°. 
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[OIV] 2D.9mu [Ke. Ill] •h.6mj [Sill] 18.7mu 




Figure 7. Predicted projected images of NGC 3918 in thi'ee infrared fine-structure lines observed by the 750 SWS, for the spindle-like model B. The upper 
panels show contour plots of the ISO-SWS apertures' transmission efficiency superimposed onto the map in the given emission line. The lower panels contain 
the projected maps of the nebula in each line after convolution with the aperture response maps. All the images are oriented such that north is up and east to 
the left. 



Table ^| provides a summary of the line fluxes predicted 
for the six infrared fine structure lines, first for the whole neb- 
ula, i.e. without the aperture corrections (column 2) and then after 
the corrections for the aperture profile, offset and orientation (col- 
umn 3). FigureQillustrates the effects of the aperture corrections 
on the final transmitted flux for three of the the six infrared fine- 
structure lines for the spindle-like model B of NGC 3918, namely 
[O IV] 25.9/im, [Ne III] 15.6/im and [S III] 18.7^m. The upper 
panel of each pair shows a map of NGC 3918 in a particular in- 
frared fine-structure emission line, oriented such that north is up 
and east is to the left; a contour map of the appropriate wavelength- 
dependent aperture transmission profile is superimposed on each 
map. The lower panel of each pair shows the map from the cor- 
responding upper panel after convolution with the appropriate 
ISO SWS aperture profile. 

It is clear from the images, and from the corrected line fluxes 
reported in Table fTsl that the SWS aperture offset is indeed the 
main reason for the SWS infrared line flux discrepancy. However, 
even after these corrections the fluxes for the predicted [Ne V] fine 
structure lines are still a factor of two to three too strong compared 
to the observed values; this may be due to the collision strengths 
bylLennon & Burke ( 199% used in the Mocassin code, which are 
much larger than those used in the Harrington code by C87. This 
has already been discussed in Section lzTI A discrepancy also re- 
m ains between th e SWS flux in the [S IV] 10.5 /im line measured 
bv lBoweJ i200l h and the corrected fluxes predicted by Mocassin, 
which are a factor of 2 to 3 stronger. C87 also predicted a flux for 
this line which was a factor of 3.6 stronger than the IRAS LRS mea- 
surements of Pottasch et alj tl986) . The origin of this discrepancy 
still remains unclear. 

It is worth noting that although the aperture offset of 14 arc- 



sec al so applied to the LWS observations made with ISO iLiu et all 
l200ll) . this would not affect the integrated flux measurements, 
due to the much larger aperture size of the LWS (~70~80 arc- 
sec FWHM). 



5 VISUALIZATION OF THE MODEL RESULTS 

A large amount of observational data has been obtained for 
NGC 3918 over the years. Some of the most recent data published 
include the optical images and the high resolution, long slit spec- 
tra presented by C99. In the same paper they also presented an 
image of NGC 3918, retrieved from the HST archives, and ob- 
tained in 1995 with the Wide Field Planetary Camera 2 (WFPC2) 
in the F555W filter and shown here in Figure [8| C99 claimed 
that the emission in this image (central wavelength = 5252 A, 
FWHM = 1223 A) is dominated by [O III] 4959-I-5007A, however 
many other nebular lines, including H/3, He I 5876 A, He II 4686A, 
[Ar IV] 4711,4740 A and [N II] 5755 A, fall within the transmis- 
sion range of this filter and , although weaker than the [O III] lines, 
can be important at particular positions in the nebula. 



5.1 Projected Images 

In order to compare the final model results for the biconical dis- 
tribution and for the spindle-like distribution, ID L visualization 
routin es, described by Morisset et al. ( 20()oh and iMonteiro et alJ 
( l200(ll . were employed. Figure |9| shows predicted maps of 
NGC 3918 in the H/3 emission line, obtained from the three mod- 
els developed in this work; namely, the biconical model (left-hand 
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Figure 9. Projected maps of NGC 3918 in H/3 for the biconical model (left-hand panels), the spindle-like model A (central panels) and for the spindle-like 
model B (right-hand panels). The top panels are displayed with linear intensity scales and the bottom panels are displayed on logarithmic scales. The spatial 
units are pixels and the intensity range is in arbitrary units. Note that a different greyscale has been used for the the spindle-like model A maps, in order to 
highlight some of the features. 



panel) and the spindle-like models A (central panel) and B (right- 
hand panel). The bottom panels of Figure|9|display the maps on a 
logarithmic scale. 

The inclination and spatio-kinematical parameters of 
NGC 3918 were derived by C99 using the spatio-kinematical 
model of ,Solf & Ulrich (.1985) . They obtained the best fit to the 
observed shape and kinematics of the inner shell by assuming an 
inclination angle S of 75° , where 9 is defined as the angle between 
the polar axis and the line of sight. The data cubes calculated in 
this work for the spindle-like models A and B were, therefore, 
rotated through the same angle to produce projected maps in 
various emission lines, including those shown in Figure |9| C87 
derived some kinematical information from a Coude spectrum 
of the [O II] AA3726, 3729 doublet. From this they derived an 
inclination angle, 9, for the biconical density distribution model, 
of 15°. Therefore the projected images for the biconical model 
(including the H/3 maps shown in Figure |3 were obtained by 
rotating the grid by the same angle. The images of NGC 3918 
published by C99 also show some low intensity small scale 
features associated with this object. In our work, however, no 
attempt was made to reproduce these features. They are discussed 
by more detail by C99, including several low-ionization structures 
located roughly along the major axis of the nebula and clearly 
visible in the logarithmic [N II] image (Figure 1 of C99, bottom 
right), where they are labeled with the letters A, A', A", B, C and 
C". The shape of the inner shell, described in Section Q is more 
clearly visible in the H ST images (e.g. Figure |8}, while the lower 
angular resolution ATT image in [N II] (Figure 1 of C99, bottom 
left) shows an almost filamentary distribution. 

From the projected maps obtained in several emission lines for 
the biconical density distribution model (the H/3 maps are shown in 
the left-hand panels of Figure|9j, it became obvious that the model 
proposed by C87, although successful in reproducing the emission 
spectrum of NGC 3918, cannot be used to actually describe the ge- 



ometric structure of this object. C87, however, were right to assume 
that this object must be optically thick in some directions and op- 
tically thin in others, in order to obtain such a complex spectrum, 
exhibiting strong emission lines from a wide range of ionization 
stages, e.g. [Ne V] and He II to [O I] and [N I]. A number of com- 
binations of thick and thin phases could possibly have been used 
to construct models to fit the emission lines and, without additional 
spatial information, it would have been virtually impossible to dis- 
tinguish between them. 

Projected maps were also obtained in several emissio n lines 
for the density distribution given by the Mellem^ il99d) model 
(spindle-like model A), discussed in Section|3| Once again, for the 
sake of space only the H/3 maps are shown in Figure |9| (central 
panel). As discussed in the previous section, the original distribu- 
tion had to be enhanced in the equatorial region in order to cor- 
rectly reproduce the observed spectrum. Although the shape of the 
inner shell of NGC 3918 is better reproduced than by the biconical 
model, the polar protrusions were still missing and the enhance- 
ment seemed to be too localised in a small equatorial torus. 

The best matches to the observed images of NGC 3918 were 
obtained for the spindle-like model B, where the density distribu- 
tion of the inner shell was reproduced using the analytical expres- 
sion given in Appendix A. Figure |9| (right-hand panels) shows the 
H/? maps obtained for this model on a linear scale (top-right) and 
on a logarithmic scale (bottom-right). The detailed shape and den- 
sity distribution given by the analytical description depends on ten 
parameters and it is possible that a number of combinations of these 
could be used to reproduce the observations. However the combi- 
nation used in this model (see Table|9j produces images (Figure|9j 
and spectroscopic results (see below) which are consistent with the 
observations. It might be possible to obtain even better fits by fur- 
ther adjusting the density distribution parameters and, possibly, also 
the central star parameters, considering that those obtained by C87 
were based on a completely different geometrical model. The cen- 
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Table 15. Predicted infrared line fluxes from NGC 3918 corrected for the 
ISO SWS beam profile and pointing position. 



Biconical Model" 


Line 


Model predictions Observed'' 




Uncorr. Corr. 



[O IV] 25.9 urn 
[Ne III] 15.6 /^m 
[Ne v] 14.3 Aim 
[Ne v] 24.3 Aim 
[S III] 18.7 Atm 
[S IV] 10.5 Aim 



352 
131 
286 
220 
48 
220 



96.1 
33.6 
69.7 
60.1 
14.3 
86.1 



89.0 
46.3 
19.6 
23.1 
8.7 
35.0 





Spindle 


like A*" 




Line 


Model predictions 


Observed'' 




Uncorr. 


Corr. 




[0 IV] 25.9 Aim 


352 


100 


89.0 


[Ne III] 15.6 Aim 


122 


21.8 


46.3 


[Ne v] 14.3 Aim 


268 


68.9 


19.6 


[Ne v] 24.3 Aim 


170 


49.4 


23.1 


[S III] 18.7 Aim 


27 


5.9 


8.7 


[S IV] 10.5 Aim 


255 


91.5 


35.0 




Spindle 


■like 




Line 


Model predictions 


Observed'* 




Uncorr. 


Corr. 




[0 IV] 25.9 Aim 


364 


104 


89.0 


[Ne III] 15.6 Aim 


126 


40.8 


46.3 


[Ne v] 14.3 Aim 


256 


58.0 


19.6 


[Ne v] 24.3 Aim 


146 


40.9 


23.1 


[S III] 18.7 Aim 


27 


10.0 


8.7 


[S IV] 10.5 Aim 


302 


114 


35.0 



The fluxes are given on a scale relative to H/3 = 100, where 

/(H/3) = 2.33 X 10-10 ergscm-2s-i. 
The fluxes are given on a scale relative to H/3 = 100, where 
/(H/3) = 2.20 X 10-1" ergs cm-2s-i. 
'^ The fluxes are given on a scale relative to H/3= 100, where 
/(H/3) = 2.56 X 10-10 ergs cm-2s-i. 



"* The fluxes are given on a scale relative to H/3 = 100, where 
/(H/3) = 2.34 X 10-10 ergs cm-2s-l 



tral star parameters, however, were based on global nebular proper- 
ties which shouldn't change as much. 



5.2 Position-Velocity diagrams and line profiles 

As menti oned above, C99 use d the spatio-kinematical empirical 
model of lSolf &Ulrici] il985h to describe the nebular expansion 
velocity, Vcxp, in order to fit their observations of position velocity 
plots of Ha and the shape of the inner shell of NGC 3918. In this 
model Vexp increases from the equatori al plane towards the polar 
axis according to the following equation ISolf & Ulrichll98^ 

Kxp(<^) = K + (Vp - K)sin^ I <?i I (4) 

where (j> is the latitude angle, defined such that (j)=Q° in the equa- 
torial plane and 0=90° in the polar direction; Vp and Vc are the 
equatorial and polar velocity components, respectively {Vp = 50 
km s-^ and 14 = 23 km s-^ in their work); 7 is a shape param- 
eter (7=12 in their work). Here the velocity field described above 



HSTF555W 





Figure 8. HST and NTT images of NGC 3918. The images are rotated so 
as to have the main jet-like feature oriented along the vertical direction. 
Left: Displayed to show the highest intensity levels by using a linear scale. 
Right: Displayed to show the faintest structures on a logarithmic scale. 
ICorradi et alll999h . 



was used for the spindle-like density distribution model A in or- 
der to produce position- velocity diagrams along slits positioned at 
P.A.=-76° and P.A.=-l-13°, in order to compare with Corradi et al.'s 
observations, shown in Figure^] Position-velocity diagrams were 
obtained for each of the lines with long slit spectra presented by 
C99, and include He II A4686, Ha and [N II] A6584. These are 
shown in Figure fTTI fmiddle two rows). 

Position-velocity maps were also produced for the biconical 
density distribution model and for spindle-like model B (Figure fTTI 
top two rows and bottom two rows, respectively) using the same 
slit positions as above. A simple radial velocity field was adopted, 
namely 



V{r) = Vi, + Vi • r 



(5) 



where r is the radial distance; Vq and Vi and constants 
(Vq = 23 km s-i and Vi = 50 km s-i in this work). Originally, the 
ISolf&Ulricili 19851) description of the velocity field was also used 
for spindle-like model B; however it was found that better matching 
results were obtained by using the simple radial velocity field de- 
scribed above, whic h was applied to obtai n the echellograms shown 
in Figures im The lSolf & UlriclJ il985h kinematical model does 
not seem to work well with the analytical density distribution used 
for model B. The disagreement is particularly evident in the major 
axis position-velocity diagrams, where the polar enhancement of 
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Figure 10. ATT long-slit spectra of NGC 3918. Images are in a logarithmic intensity scale with different cuts for the He II 6560A, Ha, and [N II] 6584A (C99) 



the velocity field creates two bright spots wh ich have not been ob- 
served. The same does not happen when the Self & Ulrich ( 1985) 
velocity field is applied to the spindle-like model A, since the two 
polar protrusions (which are created in model B by embedding an 
elliptical shell within a spheroid) are not closed at the poles (see 
Figure 0. However, from the HST image shown of the left hand 
panel in Figure |8| it seems that the polar protrusions are indeed 
closed, although not very dense, and therefore the density distribu- 
tion used in model B seems more appropriate. 

After inspection of all the position- velocity diagrams obtained 
for the three models (Figure fTTl and comparison with C99's Fig- 
ure 4, it is clear that the observations are most closely reproduced 
by the spindle-like model B (bottom two rows in Figure fTTl . All the 
synthetic position-velocity diagrams are plotted on a logarithmic 
intensity scale. It is worth noting at this point that all minor axis 
position-velocity diagrams (bottom panels in each pair of rows) 
are plotted such that the ordinate is the velocity (in [km s~^]) and 
the abscissa is the spatial variable; this has to be taken into ac- 
count when comparing the synthetic spectra to the observations in 
Figure [Tol where the abscissa and the ordinate are reversed. The 
[N II] position velocity diagrams show the largest discrepancies 
and will be discussed below. Before that, however, we will focus 
on the Ha and He II position-velocity diagrams. The projected 
maps (Figure |9j had already anticipated, and the synthetic echel- 
lograms in Figure fTTI ftop two rows) confirmed, that the biconical 
distribution is not a realistic representation of NGC 3918. How- 
ever, as already mentioned, the imaging information available for 
NGC 3918 was very limited at the time when this model was de- 
vised by CSV. Model A (middle two rows of the position-velocity 
diagrams in Figure fTTl fails to reproduce the polar protrusions in 
the observed long slit spectra (Figure fTol . as was also clear from the 



projected maps obtained from this model (Figure|9| central panels). 
Moreover, the thick waist in model A, which had to be artificially 
created in order to obtain the observed line strengths for the lower 
ionization species, is not observed in NGC 3918. 

None of the models presented in this work succeeded in repro- 
ducing the [N II] A6584 echellograms shown in Figure [Tol There 
are two main reasons for the observed discrepancies: firstly, N"*" 
is abundant in the FLIERs of NGC 3918 (C99), which are not ac- 
counted for in the models presented here; secondly, N'^ also seems 
to exist in the inner regions of the shell; the NTT image in Fig- 
ure 1 of C99 shows [N II] emission to be present in almost a fila- 
mentary form. All the models presented in this work have smooth 
density distributions and do not include any filaments or knots. In 
spindle-like model A, N'^ is mainly concentrated in the toroidal 
region behind the thick equatorial waist, which is screened from 
direct starlight, whereas in model B, N"*" is mostly concentrated 
in small polar regions beyond the elliptical shell which creates 
the protrusion, and they are also screened from direct starlight by 
this shell and are further away from the star. In fact, it is quite 
hard to have N"*" survive in the very hard stellar radiation field 
(Teff = 140,000 K) in the inner shell, unless some sort of density 
enhancement is postulated. This shows that, although spindle-like 
model B appears to be quite successful in reproducing the main 
spectroscopic and kinematic characteristics of NGC 3918, there is 
still scope for improvement. 

A three-dimensional reconstruction technique for studying the 
morphology, physical condition, ionization and spatial structure of 
planetary nebu l ae has been introduced by Sabbadin et al. (200(|), 
ISabbadin et alj i200(j) and iRagazzoni et ah (.2001.) . In this app- 
proach, the relative density distribution of an emitting region along 
the cross-section of nebula covered by the slit can be obtained 
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Figure 11. Position-velocity diagrams in He II 4686 A (left-hand panels). Ha (central panels) and [N II] 6584 A (right-hand panels) for tlie tlie biconical 
model (top two rows) and spindle-like models A (middle two rows) and B (bottom two rows) of NGC 3918. The top panels in each pair of rows show the 
echellograms obtained with slits oriented at PA = -76° while the bottom diagrams show echellograms obtained in the same lines, but with slits positioned at 
PA = -1-13°. 
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from the radial velocity, the FWHM and the intensity profile. F. 
Sabaddin and other members of the Padova group have recently 
obtained observations for several other planetary nebulae, includ- 
ing NGC 3918, consisting of long slit echellograms at various po- 
sition angles on the nebulae, and are planning to carry out detailed 
studies of them by applying their three-dimensional reconstruction 
technique to them. This means that in the near future a data cube 
containing a three-dimensional density distribution obtained using 
the method above could be available for NGC 3918. It would be 
very interesting then to map such a data cube onto a Mocassin grid 
and construct a new photoionization model for this object, since as 
discussed earlier, although spindle-like model B produces satisfac- 
tory results, there are still some discrepancies remaining with the 
observational data. 



6 CONCLUSIONS 

In this paper three photoionization models were constructed for the 
planetary nebula NGC 3918, a biconical model and two spindle- 
like models. T he first spindle-l ike model (A) used the density dis- 
tribution of the lMellemaHl99(i) model already applied to this object 
by C99. The second spindle-like model (B) instead used an analyt- 
ical expression for the density distribution, which aimed to mimic, 
by means of an ellipsoid embedded in a sphere, the shape of the 
inner shell of NGC 3918. The aim was to find a model which could 
not only reproduce the main spectroscopic features observed for 
this object, but also the spatio-kinematical properties of the nebula 
which recent observations (C99) have uncovered. 

The integrated emission line spectra obtained from the three 
different models were all in fair agreement with the observa- 
tions. Discrepancies exist between the observed line fluxes of the 
C IV AA 1548, 1550 and Mg IlAA2796,2803 resonance doublets, 
which can largely be e xplained by dust absorption in the nebula 
iHarrington et alj[l988h . Large discrepancies were also found be- 
tween the models' predictions for the fluxes of the infrared fine- 
structure lines and the measurements of the ISO SWS spectra. The 
main cause for this discrepancy is a pointing error which affected 
ISO observations of NGC 3918, causing an offset of approximately 
14 arcsec from the centre of the nebula. The corrected predicted line 
fluxes, obtained by convolving the ISO SWS aperture profiles with 
projected nebular maps in the relevant emission lines, confirmed 
that most of the discrepancy between the observed and model fluxes 
could be attributed to the pointing error, although there is still a fac- 
tor of three discrepancy for the [S IV] 10.5 /im line. 

A diffuse radiation field consistency test was also carried out 
in this work, which showed that the interaction of the diffuse radi- 
ation fields from two adjacent regions of different densities is not 
negligible, even in the relatively uncomplicated case of the biconi- 
cal density distibution used by C87. We found that the low ioniza- 
tion species in the optically thick cones were particularly effected 
by the diffuse radiation field coming from the optically thin sector. 

Although the volume-integrated emission line spectra ob- 
tained by the three models of NGC 3918 were in agreement with 
each other, the projected maps and the synthetic long-slit spectra 
obtained in several emission lines were, however, very different 
from one model to the other Spindle-like model B produced the 
best fits to the observations, although some discrepancies still ex- 
ist, particularly in the [N II] maps and long slit spectra, as discussed 
in Section l5!2l 

Confirming the conclusions of iMonteiro et alj i200(]h . from 



their three-dimensional modelling of NGC 3132, this work has in- 
dicated that a detailed model of a nebula cannot be verified just by 
comparison of the observed integrated spectrum with model pre- 
dictions. In fact, in the case of NGC 3918 approximately the same 
spectrum can be obtained with a number of different geometries 
and density distributions. For this reason, three-dimensional mod- 
els are necessary in order to allow a spatio-kinematical analysis to 
be carried out by comparing predicted images and position-velocity 
diagrams in several lines to available observational data. Mocassin 
provides all the tools needed for such simulations and for the visu- 
alization of the final results. 
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Appendix A: Analytical description of spindle-like model B 

We define the following cut functions: 



, , 1 { X — Xo\ 

Coi(a;,a;o, At) — 0.5 H tan — 

TT V Aa; / 

Ciofa;, xo, Ax) — 0.5 — — tan~^ ( ^ . ^° 1 

TV \ Ax J 



(6) 



The functions above can be used to describe a smooth transi- 
tion from to 1 and from 1 to zero, respectively, about the value 
xo, with a typical size for the transition of A x. 

The equatorial spherical shell can then be described as fol- 
lows: 

equat{r, 9) = Coi{9, 9i, A6)-Coi{r, r^q, Ar"q)-Cio{r, r^q, Ar2q^){T) 

where r is the radial distance and 9 the polar angle (9 — 90° 
at the equator and ^ = 0° at the two poles). Values for 9i, A9, Veq, 
Ari"^ and Ar°f are given in Table|51 

We then define r^, the radius of the ellipsoid, as 
Te — \/ x'^ + y'^ + {e ■ 2)2, where e is the eccentricity (given 
in Table|9), then the polar protrusions can be created using the fol- 



fluxes of Ne''+, namely [Ne v] 14.3 /im and [Ne v] 24.3 fim, which 
disagree by factors of almost 10 and 5, respectively with C87's 
predictions and indeed with the observations. The discrepancy be- 
tween the models is due to differences in the collision strengths 
used for these lines. The collision stren gths for these transition s 
used by Mocassin are the calculations bv lLenno n ^^u^^ ^9^|), 
based on the R-matrix method described bv iBeMngtor^^l] 
1 1987). These collision strengths are much larger than t he rates 
gi ven by Mendoza ( 1983.) based on the results of Saraph et alJ 
al97d) . which were used by C87 in their mo del. In 1987, how- 
ever, the collision strengths o f lAggarwall iT98l were already avail- 
abl e for [Ne v], and these were also much higher than those given 
by iMendozj il983h . since they took into account the large res- 
onances existing at energies just above the threshold for excita- 
tion of the fine structure lines. C87 decided against using those 
results as they were not in agreement with the observed ratios of 
[Ne v] 3426 A/ 14.3 /im. However, van Hoof etal. |2000) later 
argued that the collision strength s calculate d by Len non & Burk3 
1 1991), which were similar to the lAggarwaL (,1983.) values, are re- 
liable, and that the discrepancies found by C87 can be explained 
by the inaccuracy of the [Ne v] 3426 A flux that they adopted. 
Moreover, Berrington |2001) in his contribution to the November 
2000 Lexington workshop, discussed near-threshold resonances in 
the calculation for collision strengths for fine-structure transitions, 
taking [Ne V] as a case study. He compared the effective collision 
strengths obtained using different methods, showing that they all 
largely agree with the Lennon & Burke 1 1991) results. During the 
same Lexington meeting, a new R-matrix calcu lation of [Ne V] was 
independently published bv lGriffi n & Badnelll iJoOQ) who also ob- 
tained similar results . In the cur rent work the decision to use the 
collision strengths of lLennon & B urke ( 1991), although producing 
results which are in disagreement with C87's observations, was fi- 
nally taken, since the benchmarking carried out for the Lexington 
optically thick planetary nebula, which has a central star with ef- 
fective temperature, Te// — 150, 000 K (similar to the central star 
of NGC 3918), clearly showed that most photoionization codes 

all but Sutherland's Mappings) use the new collision strengths 

lErcolano et all2002h . 

The likely cause of the large discrepancy between the pre- 
dicted and observed [Ne V] fine-structure line fluxes is addressed 
in Section|4]and in Appendix C below. 



Appendix C: The ISO pointing error 



lowing expression 

pole{r, 9) = Cio{9, 9i, A9)-Coi{re, rpo, Arp")-Cio(re, rpo, Ar"" 

where rpo, Ar^ and Ar^" are also given in Tabled Finally, 
the density distribution at any point in the grid can be obtained by 
combining EquationsQand|8| as follows 

NH{r,9) = nix{equator{r,9) + pole{r,9)). (9) 



Appendix B: The [Ne v] flne-structure lines 

The largest discrepancies between this work and that of C87 
were obtained for the predictions of the infrared fine structure line 



Most of the infr ared line fluxe s measurements reported in Ta- 
ble|2|are taken from lBoweJ i200ll) (reference d in Table|2}, mea- 
sured from spectra taken with the ISO Short- Wavelength Spec- 
trometer (SWS) between 1996 and 1998. The fluxes measured by 
iBoweri C2001') for the infrared lines of NGC 3918 are systemati- 
')(8)cally sm aller, by approx imately a factor of two, than those mea- 
sured by IPottasch et alj ll9S(h. who used IRAS Low-Resolution 
Spectrometer (LRS) observations and are smaller by factors of 3 to 
6 for most lines, apart from the [Ne V] fine structure lines, which 
are between 10-15 times smaller, compared to both Mocassin and 
C87 predictions. The reason for the [Ne V] fine structure line dis- 
crepancy between the Mocassin and C87 predictions has already 
been discussed above. This ISO SWS discrepancy is also evident 
in the spindle-like models of NGC 3918 (see Section ITTl and is 
believed to be due to the combined effects of incorrect object co- 
ordinates used for the ISO observations and the small size of the 
ISO SWS aperture. The ISO LWS fluxes for the [N III] 57.3 ^im 



and [O III] 51.8 /im and 88.4 fim lines are in good agreement 
with the model predictions. Although the same pointing error of 
14 arcsec occured, the much larger LWS aperture of 70-80 arc- 
sec minimised the effects on the observed fluxes. Corrections to 
the predicted ISO SWS line fluxes are calculated and discussed in 
detail in Section|4] 

The position given by the Strasbourg-CDS Simbad database 
for NGC 3918, in the J2000 FK5 reference frame, is 
q(2000) = 11 50 17.2 and a(2000) = -57 10 53, while the posi- 
tion used for the IS O 5^5 and LW S observations of NGC 3918 
was that measured bv lMilneH l97^. namely a(1950) = 11 47 50.1 
and (5(1950) = -56 54 10, equivalent to a(2000) = 11 50 18.9 and 
(5(2000) = -57 10 51. This resulted in the SWS apertures being off- 
set from the nebular centre by 1.8 arcsec north in declination and 
by 13.6 arcsec (= 1.68 sec) east in right ascension, giving a total 
offset of 13.8 arcsec from the centre of the nebula. The effects of 
this SWS aperture offset, together with the particular spacecraft roll 
angle at the time of the the observations (retrieved from the ISO 
archives) and the SWS subspectra beam profiles, are sufficient to 
largely explain the discrepancies initially found between the model 
predictions and the observations. The wavelength dependent aper- 
ture sizes of the ISO SWS range from 14"x20" for 2.38/im^ 
A ^ 12.0^m, 14"x27" for 12.0^im A 27.5^im, 20"x27" for 
21.5um A < 29.0^m , to 20" x 33" for 29.0/im ^ A 45.2/im 
<deGraaw&etalJll99d) . Five of the six infrared fine-structure 
lines listed in Table 1 151 were observed by the ISO SWS through a 
14"x27" aperture, while the [S IV] A10.5/im line was observed 
through a 14" x20" aperture. As mentioned above, the aperture di- 
mensions are not the only factor to be taken into account, there is 
also the transmission efficiency profi le across the spatial a perture, 
which, in some cases, is far from flat Leech & et all2002h . 



3-D photoionization modelling of NGC 3918 



